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Abstract 

We study both Bose and Fermi gases at finite temperature and density in an approxi- 
mation that sums an infinite number of many body processes that are reducible to 2-body 
scatterings. This is done for arbitrary negative scattering length, which interpolates be- 
tween the ideal and unitary gas limits. In the unitary limit, we compute the first four 
virial coefficients within our approximation. The second virial coefficient is exact, and 
we extend the previously known result for fermions to bosons, and also for both bosons 
and fermions for the upper branch on the other side of unitarity (infinitely large positive 
scattering length). Assuming bosons can exist in a meta-stable state before undergoing me- 
chanical collapse, we map out the critical temperatures for strongly coupled Bose-Einstein 
condensation as a function of scattering length. 
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I. INTRODUCTION 



The growing amount of increasingly accurate data from experiments on cold atoms 
[l-6] poses particularly interesting challenges for theorists to develop new methods. 
This is especially true for experiments where the scattering length can be tuned 
to vary anywhere between ±00 using Feshbach resonances. Monte-Carlo methods 
have been develo ped sufficiently that excellen t ag reement with experiments has now 
been achieved [l-ioj]. Recent reviews are 11, 12]. Nevertheless, the development of 



new analytical methods, though approximate, continues to be a worthwhile pursuit 
because they can afford new insights into the underlying many-body physics. 

One such method has been developed by one of us, and is based entirely on the 
zero temperature S-matrix 131 . It is reminiscent, in fact was modeled after, the 

n 

thermodynamical Bethe ansatz equations of Yang and Yang [14[. The ingredients 
are the same: the occupation numbers are parametrized in an ideal gas form, but 
with the one-particle energy replaced by a pseudo-energy e(k). The latter satisfies 
an integral equation with a kernel based on the logarithm of the 2-body S-matrix 
at zero temperature, and there is a simple expression for the free energy at finite 
temperature and density. Whereas for integrable theories in 1 spatial dimension the 
thermodynamic Bethe ansatz is exact because of the factorizability of the many- 
body S-matrix, the formalism in 13J is certainly an approximation. Nevertheless 
it has certain desirable features, such as the fact that the 2-body S-matrix can be 
calculated exactly in non-relativistic theories, and has been demonstrated to give 
reasonably good results in some regimes. For instance, it was applied to the so-called 
unitary limit in 3 dimensions where the scattering length diverges, and the S-matrix 
becomes simply —1, and reasonable results were obtained for the critical temperature 
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16| . The ratio of the viscosity to entropy density was also ca 



method [n| and agrees well with the most recent experiments 18 



culated using this 
. Thus, although 
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the method cannot really compete with numerical methods such as Monte-Carlo, it 
can be justified as an exploratory tool for regimes that have not been extensively 
studied. 

This paper is mainly concerned with using the method to study the critical prop- 
erties of Bose and Fermi gases in 3 dimensions as a function of scattering length, in- 
cluding the vicinity of the unitary limit where it diverges. For 2 component fermions 
this is the familiar BEC/BCS cross-over. For negative scattering length the interac- 
tions are attractive and there is a phase transition to a strongly coupled version of 
superconductivity. For positive scattering length the fermions have a bound state, 
i.e. the 'atoms' form 'molecules', which can subsequently Bose-Einstein condense. 
This fermionic case has been already extensively studied and we have nothing novel 
to report here. 

On the other hand the bosonic case has been much less studied theoretically and 
is just beginning to be explored experimentally, and this is the main subject of this 
article. The spectrum is analogous to the fermionic case: for negative scattering 
length there is no bound state, whereas for positive scattering length molecules can 
form via 3-body processes. We thus will restrict our study to the case of negative 
scattering length. This case has perhaps not been studied very much theoretically 
because it is believed that the attractive interactions lead to a mechanical instability, 
i.e. the gas collapses. However it remains possible that this state could exist as a 



meta-stable one 19]. 



On the other side of unitarity, in our analysis we would need to incorporate the 
molecules, with their own pseudo-energy etc, and this is beyond the scope of this 
work. However for the so-called 'upper branch', the molecules are assumed to be 



absent, and this situation has been studied experimentally 20M22| and theoretically 
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26|. This motivated us to present new results on the virial expansion for both 



fermions and bosons on this upper branch. 



Our results are presented as follows. In the following section we review the S- 
matrix and renormalization group for the models and present our conventions for 
the coupling and its relation to the scattering length, which are the usual ones. In 
section III we review our method, describing in a precise way what we are neglecting 
in the approximation, and how in principle to calculate the corrections. The virial 
expansion is studied in section IV, where we reproduce the known second virial 
coefficent for fermions on the BCS side, but also include new results for bosons and 
for the upper branch. Here we also calculate the third and fourth virial coefficients 
in our approximation in order to compare them with more accurate calculations and 
experiments. In section V we study the extension of the Bose-Einstein condensation 
of the ideal gas to the full range of negative scattering length, thereby mapping out 
how T c depends on the scattering length. In section IV we revisit the fermionic case, 
extending the results in [igJ to arbitrary negative scattering length. 



II. CONVENTIONS, S-MATRIX, SCATTERING LENGTH 

The bosonic model we consider is defined by the action for a complex scalar field 

0. 

d^dtL^^- 1 ^- 9 -^^) (i) 

where positive g corresponds to repulsive interactions. By the Galilean invariance, 
the two-body S-matrix depends only on the difference of the incoming momentum 
of the two particles k, k': 
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Unitarity of the S-matrix amounts to S*S = 1. 

The momentum space integrals for the higher loop corrections are divergent and an 
upper cut-off A must be introduced. In the above expression, g R is the renormalized 



coupling: 

1 1 rnA 

Denning g = 'g/A, where g is dimensionless, and requiring gn to be independent of 
A gives the beta-function: 

where I = —log A is the logarithm of a length scale. The above beta function is 
exact since it was calculated from the exact S-matrix. One thus sees that the theory 
possesses a fixed point at the negative coupling g* = — 27r 2 /mA where it becomes 
scale invariant. 

We turn now to the scattering length a s . It can be defined as j^ : 

lim — = -a s (5) 

|k-k'|->o |k — k'| 

where 5 = — i log S. From the above expression for the S-matrix, one finds 5 = 
— 2 arctan(mg/j|k — k'l/Svr), which gives 

mg R f 

as = ^r (6) 

One sees that scattering length diverges at precisely the fixed point g = g*. Heuristi- 
cally, the loss of this length scale implies universal properties of the free energy since 
it can only depend on the chemical potential and temperature. Note the S-matrix 
becomes S — — 1. The scattering length a s — > ±oo, depending on from which side 
g* is approached. When g = g~, i.e. just less than g*, then a s — > oo, whereas when 
g = gf, a s — > — oo. For reasons described above, this paper will mainly only consider 
negative scattering length — oo < a s < 0. 

For fermions we consider the two-component model defined by the action: 
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(7) 



With this convention for g, the S-matrix is the same as for bosons, eq. ([2]), as is 
the beta function and scattering length eqs. ( Hf6|) . For fermions, negative scattering 
length corresponds to the BCS side of the BCS/BEC crossover; we will thus only be 
working on the BCS side. 

For positive scattering length, the S-matrix has a pole signifying a bound state, 
or "molecule". In order to have a smooth crossover across the unitary limit, this 
bound state must be incorporated into the thermodynamics, and this is beyond the 
scope of this paper. Thus we will be primarily studying negative scattering length 
where there are no molecules. The system on the other side of unitarity where 
molecules are ignored is usually referred to as the "upper branch" . In certain regions 
of density/temperature, the upper branch can in fact be metastable, and has been 
realized in experiments. We will thus present a few results on the virial expansion 
for the upper branch for both bosons and fermions. 



III. S-MATRIX BASED FORMALISM FOR THE THERMODYNAMICS 



The method developed in 13] is based entirely on the S-matrix. Contributions to 
the free energy density J 7 have a diagrammatic description. Vertices with 2N legs 
represent the logarithm of the S-matrix for N — > N particle scattering. Diagrams 
that contribute to J 7 are closed diagrams with vertices linked by occupation numbers 
/o(k) = \/{eP^~^ - s) = z/(e^^ k - s z), where the temperature T = fi is the 
chemical potential, z = e M / T , cu^ = k 2 /2m is a 1-particle energy, and s = +1,-1 cor- 
responds to bosons, fermions respectively. Vertices exist for any N. These diagrams 
are not to be confused with finite temperature Feynman diagrams in the Matsub- 
ara formalism. Here, vertices represent the S-matrix to all orders in perturbation 
theory at zero temperature, i.e. each vertex already represents an infinite number 
of zero temperature Feynman diagrams; the finite temperature dependence comes 
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about mainly from the occupation numbers on the legs. This method thus appears 
inherently different from the t-matrix approach for instance. The low order diagrams 
are shown in Figure [TJ Explicit expressions for some of these diagrams will be given 
in the next section. 




(a) (b) (c) (d) 



FIG. 1: Diagrams contributing to the free energy density, where J-q is for the ideal gas. 

Two-body processes are by definition diagrams built only out of 4-vertices, such 
as those in Figures l(a,b,c). There are an infinite number of diagrams contributing 
to just the two-body processes. An infinite subset of two-body diagrams are so-called 
foam diagrams, such as in Figure l(a,b), but with an arbitrary number of bubbles. 
These can be resummed from a variational principal based on the diagram Figure 
la. The result is the following. 

It is convenient to express the free energy density J 7 and density n in terms of 
scaling functions c, q of the dimensionless variables x = fi/T and a = \t/cl s -, where 
Xt = y^27i/mT the thermal wavelength, as follows: 

n\% = q(x,a) (8) 
JTA| = _((5/2)Tc(x,a) (9) 

where £ is Riemann's zeta function. With the above normalizations, for the ideal 
gas at zero chemical potential, q = C(3/2) and c — 1. The two scaling functions c 
and q are of course related since n = —dJ^/dfi, which leads to q = ((5/2)d x c. 

For fermions, is convenient to define the Fermi surface wavevector kp = (S^n) 1 ^ 3 , 
where n is the 2-component density, and Tp = k F /2m. In terms of the single 
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component scaling function q: 

2/3 



T ( 4 



(tt^) ' r- = -(6A)- 1/3 (io) 
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The definitions leading to the above formulas make sense also for bosons; for instance, 
the BEC transition of the ideal Bose gas occurs at T C /T F = (A/3 v ^r((3/2)) 2 / 3 = 
0.4361. 

In our formalism, the filling fractions, or occupation numbers, are parameterized 
in terms of a pseudo-energy e(k) in an ideal gas form: 

f d 3 k 1 

n J (2tt)3 e^)/ T - 1 1 ' 

The pseudo-energy can be thought of as a 1-particle energy in the presence of all 
the other (interacting) particles in the gas. The consistent summation of many body 
processes that involve only 2-body scattering described above leads to an integral 
equation for the pseudo-energy e(k), analogous to the Yang- Yang [l^] integral equa- 
tion. It is convenient to define the quantity: 

y(k) = e -(eM-«k+/*)/T ( 12 ) 

Then y satisfies the integral equation 

If d 3 k> yjk')- 1 
y^) = l + Tj j2^ Gs{k - k) e<«)/T-s (13) 

where refers to boson/fermion. 

The kernel G is related to the logarithm of the 2-body S-matrix of the last section: 

G+(k,k') = arctanf^lk-k'lV G_ = G+/2 (14) 

m\K — k'| V 87r / 

In the unitary limit a s — > Too, the kernel takes the simple form: 

G + (k,kO = ± * n \, (15) 
m k — k 



where the positive sign corresponds to infinite negative scattering length. 

By rotational invariance, y is a function of |k| 2 . It is convenient to define the di- 
mensionless variable k = \/|k| 2 /2mT. The angular integrals in the integral equation 
( jl~3l) can be performed analytically (Appendix A). The result is the following: 



oo 

z a 



a 2 /Tr+(k + k') 2 

it J Q e k ' 2 — szy(k') \ 2k^pii & \_a 2 /n + (k — k') 2 



V(k) = 1 + - I dk'k' -p 7 ^- 7 = log 



7T 



j + 1 J arctan ^ — (k + k') ) (16) 
k' 



1 arctan ( (k — k') 

k J V a 



The scaling functions then have the form: 



* r**M^ (17) 

Vtt J e k - sy(k)z 



and 



dkk 2 ^—s\og(l — szy(k)e~ 



k 2 



1 z(y(k) - 1) 



"V^C(5/2)y V V yw J 2e k2 -szy(k), 

The ideal, free gas limit corresponds to y = 1 where q = sLi3/2(s2;) and c = 
sLi 5 / 2 (sz)/((5/2), where Li is the polylogarithm. The BEC critical point of the 
ideal gas occurs at fi — 0, i.e. q = ((3/2). 

Since the fermion model has two components, in equations [HJ q, c — > 2q, 2c with 
eqs.( [T7]fT8|) still valid. In other words, henceforth, q, c will refer to one of the two 
components. 



We will need the entropy per particle. The expression in 16] must be generalized 
for q^O, given the extra T dependence in J 7 through a. One has 

<9.F(/i, T) dT dx dT da dT^i, T, x, a) 



dT dx dT da dT dT 



(19) 
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where in the last derivative only the explicit T dependence is considered. This gives 



: C ( 2 ) X t 



-c — xc 



-ac 



2 2 

where c = dc/da and c' = dc/dx . The entropy per particle then takes the form: 



s/n = C 



q 



5 , 1 . 

-c — xc ac 

2 2 



(20) 



IV. VIRIAL EXPANSION 



The virial expansion is conventionally defined as a series expansion of J 7 , or the 
density, in powers of z: 



n=l 
oo 



(21) 



n=l 



where the second relation follows from n = —dF/dfi. In the free theory, the series 
expansion of sLi 5/2 (s2;) gives bf\ bf\bf\ bf\ ... = 1, , ^, .... Expanding the 
occupation number in powers of z, fo = ze~ k2 + z 2 se~ 2k2 + each diagram in 
Figured] can be expanded in powers of z. Since each internal leg corresponds to an 
fo, a diagram with m internal legs contributes all b n with n > m. Thus, the exact b% 
comes from Figure la. The contributions to 63 come only from Figures la and Id, 
and to 64 from Figures la, lb, lc, Id, and an additional diagram not shown which 
is "primitive" as in Figure Id, but with a vertex with 8 legs and 4 loops. 

In the two-body approximation captured by the integral equation of the last sec- 
tion, only Figures la and lb contribute to the first 4 virial coefficients. Let J 7 ^' 6 ) 
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denote the contributions to J 7 from these two diagrams. They are given by 13 1 
^ - 1/ /o(k)^(k,k')/o(k') 

T<t) = W /^|^!^s/o( k .) G -( k ..^)/»( k ^ 2G -( k - k 3)/o( k 3) (22) 
Expanding the free occupation numbers /o in powers of z, the contributions to 0i_4 
from diagrams Figure la, lb, denoted b^_ A and bf_ 4 respectively, are the following: 
XI f d 3 k d 3 k' 



h (a) _ At f d k d k /T „ /t q , f] 

b2 ~ 2Tj (27r)3(2vr)3 e 6 Gs(k ' kj 

^ = # / (SS e ^ /Te ^ /TGs(k> k0 (e "" /T + e "" /T) (23) 

4 ~ 2T J (27r)3(2vr)3 e 6 +C +C > 



whereas involves two kernels G: 
»?' = § / ^ IT ^' T ^ IT G^M)G M (24) 

A. Infinite negative scattering length 

In the unitary limit, the kernels have the simple form eq. (I15j) . and the virial 
coefficients are pure numbers by the scale invariance. The above integrals for 6 2 ,3 are 
easily performed analytically by making the change of variables ki = k — k', k 2 = 
k + ak', where a is chosen to cancel the cross term in the exponential; the result then 



factorizes into two 


gaussian integrals. 


Chosing 


a = n'/n, one can show: 




A* r 


d 3 k d 3 k' 




4vr 2 1 


(25) 


2T J 


(2tt)3 (2vr)3 6 




mk-k' nn'y/n + n' 


The result for o 2)3 


— 2 ,3 + °2,3 lb 


then: 










, 3V2 
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h = 


8^ 
27 


(fermions) 


(26) 




h = 8 ' 


b 3 = 


19^3 
27 

11 


(bosons) 





Finally let us turn to 64 which has contributions from 64 and Figures la and 
lb. The contribution b± involves only one kernel and amounts to a sum of terms 
involving the integral eq. (]25|) . giving 64 = 11/24 and 11/12 for fermions and 
bosons respectively. The integral for can be simplified by noting that the k x and 
k 3 integrals are identical. After rescaling the k's, for fermions one obtains 

6f = ^/A 2 e-i(/ rf 3 kl _£^_y (27) 

Shifting ki — 7- ki + k2, and performing the angular integrals one finds that the ki 
integral in parentheses is 7r 3//2 Erf(/c2)/^2- Putting this all together, one obtains for 
K = bf+b^+bf: 

41 f°° 

6 4 = — -2v^F/ dke~ 2k2 Eri(k) 2 = -0.0535 (fermions) (28) 

91 f°° 2 

b 4 = — + 8v^F / dke~ 2k Erf(fc) 2 = 2.870 (bosons) 
96 Jo 



It is interesting to compare with the known results for 
agrees with the exact result calculated by other methods 



ermions. As expected 62 



281 ] . whereas 6 3j4 are not 



since some 3 and 4 body physics has been neglected. The calculation in the work 



29 



30| includes 3 body physics and obtains 63 = —0.29, in very good agreement 
with experiments compared to our result 63 = —0.5132. Although the sign is 
correct, this does indicate the importance of the 3 body corrections. The coefficient 
64 has been extracted from experiments 4] with the result 64 = 0.065, compared to 
our —0.0535; here we obtain the correct trend that the absolute value of the 6's is 
decreasing, however the sign is incorrect. Note the large value of 64 for bosons verses 
fermions. 
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B. The upper branch: Infinite positive scattering length 



Let us present results for the upper branch on the other side of unitarity. Here, 
the kernel changes sign eq. (IT5"|) . For fermions this leads to 



5V2 



10V3 
27 ' 



47 f°° 

b A = 2v^F / dke~ 2k ' 'Erf (k) 2 = -0. 

96 Jo 



9702 



and for bosons: 



7^2 



17\/3 
27 ' 



or /"OO 

b 4 = -— + 8v^ / dke~ 2k2 Erf (k) 2 = 1.034 
96 Jo 



(fermions) 



(bosons) 



(29) 



(30) 



V. ANALYSIS OF BOSONS: STRONGLY COUPLED BEC 



As explained in the Introduction, in the case of a negative scattering length, we 
assume the bosonic gas can exist as a metastable state, i.e. stable against mechanical 
collapse, and has a BEC transition smoothly connected to that of the ideal Bose 
gas. Recall that for the ideal gas, the BEC transition occurs at x = 0, leading to 
q c = nAy = £(3/2). We wish to explore q c as a function of a = \r/a s - 

As for the ideal gas, the criterion for BEC is clearly the condition on the pseudo- 
energy e(k = 0) = 0, since this implies the occupation number / diverges at k = 0. 
We solved the integral equation (116j) iteratively until the solution converged. A plot of 
e(k = 0) as a function of x leads to the identification of the critical value x c (a), which 
in turn gives q c {a). We plot our results for the critical values of x c and q c = n\^ c in 
Figures El [3j The equivalent results in terms of T c /Tp as a function of l/kpa s are 
plotted in Figure HI In all these figures one sees that as a s — > 0, one recovers the 
ideal gas BEC results x c = 0,q c = C(3/2) = 2.61, and T C /T F = (4/3C(3/2))v^F) 2 / 3 = 
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0.436. In the opposite limit a s — > oo the result is also consistent with the unitarity 



limit results in 



16|. 
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FIG. 2: x c (a) for different values of a = Ar/a s for bosons. As expected, as a — > — oo, 
x c -> 0. 




FIG. 3: Critical temperature (or density) as a function of Xt/cl s for bosons. 
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FIG. 4: Critical temperature in units of Tp as a function of \jkpa s for bosons. 

VI. ANALYSIS OF FERMIONS 

For fermions with negative scattering length, the attractive force leads to Cooper 
pairing and thus a strongly coupled superconductivity that connects smoothly with 
the BCS theory at very small scattering length. In this section we attempt to un- 
derstand this phase transition within the present formalism. In this case however 
the criterion for the transition is not as obvious as for the BEC transition of the 
last section. We will pursue the hypothesis that a transition in the behavior of the 
entropy per particle can be used as a signature of the phase transition, as was done 



m 



16]. 



The entropy per particle as a function of TjTp and /i/T are shown in Figures [51 
[6l and Figure [7] relates the chemical potential to TjTp. Figure [7] indicates that the 
chemical potential \x increases with \a\, as expected for stronger interactions. 

As explained above, we use the entropy per particle plots to identify the critical 
point. One sees from Figure [5] that the entropy begins to increase as a function 
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FIG. 6: Entropy per particle as a function of x = n/T. 

of temperature when the temperature is low enough, which is interpreted as being 
unphysical; it is actually ill defined below this temperature, i.e. there is no longer 
a solution y to the integral equation. This change in behavior is also seen in Figure 
[61 For instance, in the unitary limit a = 0~, the critical point is at x c ~ 11. 
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FIG. 7: \ijtF as a function of T/Tp. 



which corresponds to T/T c « 0.08. In Figure |8] we plot our results for T c /Tp as a 
function of 1/ kpg^. Our results are in rough agreement with the most recent Monte- 



Carlo results 17 



101 ] . In the range —0.5 < l/kpa s < the Monte-Carlo results are 



0.08 < T c /Tp < 0.15, whereas we obtain 0.07 < T c jTp < 0.08, which is on the low 
side. 



VII. CONCLUSIONS 

We have applied the S-matrix based formalism developed in 13l to bosons and 

n 

fermions with arbitrary negative scattering length, extending the results in [16( at 
the unitary limit. We explored the virial expansion up to the fourth order within 
our approximation, confirming that 62 for fermions is exactly correct, i.e. agrees 
with known results computed by different methods. We extended these calculations 
to bosons and both bosons and fermions to the other side of unitarity, that is, the 
upper branch. Our value for 63 for fermions on the other hand, when compared with 
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0,085- 




1/k F a 

FIG. 8: T C /T F function of l/kpa s for fermions. 

more accurate calculations and experiments, indicates the importance of the need to 
incorporate 3-body processes, and we hope to carry this out in the future based on 
the prescription described in jl^ . 

We applied the method as an exploratory tool for less-studied regimes. In par- 
ticular, for bosons, assuming they can be meta-stable against mechanical collapse, 
we calculated T c for a strongly coupled Bose-Einstein condensation as a function of 
(negative) scattering length. We hope that our results in this case may be compared 
with other methods, or better, experiments, in the near future. 
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IX. APPENDIX A 



For simplicity consider the bosonic case s = +1. By using |k — k'| 



y/k 2 + {k') 2 - 2kk' cos# , where |k| = k, and integrating in spherical coordinates, 
one has: 



vW = i - ^ 



rriTT 



dp p 2 



de 



sin 9 



v - zy(k') J y/k 2 + p 2 -2kpcos9 



• arctan \^\Jk 2 + p 2 -2kp cos 8 



1 + 



4/3 



rrnr 



dx 



oo z n-1 



arctan ( m ® R y/a — bx 
V 8vr 



(31) 



where p = |k'|, a = k 2 + p 2 , and b = 2kp. In the third line above, the function 
y(k) is considered to be spherically symmetric (y(k) = y{k)). 
Evaluating the second integral 

1 dx + ( m 9 R I TT\ 

— ^=^= arctan y a — bx 

i Va - bx \ 87r / 

— —(k + p) arctan ( — —(k + p)) + (k — p) arctan ( — —(k — p) 
kp L V 87r / V 87r 



47T 

H log 

mg R 



one obtains: 



= 1 + / dpp 2 7T 7 < log 2 

m7r J e fip /2m _ zy ^ k 1 mgR / x + ^ _ ^ 2 

— (fc + p) arctan ( s f R (k + p) J + (fc — p) arctan ( ® R {k — p) 

V 87T / V 87T 



f _8tt_ 



+ (A; + p) s 
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Rescaling k — > \/2mTk, p — > y/2mTp, and introducing a = Xt/cl s , one finds that the 
integral equation of y(n) can be written as in eq. (|T6|) . 

The eq. f lT6|) is valid for any value of g^. As a consistency check, one can verify 
that in the unitary limit gn — > ± oo =>- a — > O^ 1 , one obtains 



lzp 4 / oIk' — 



e K ' — zy(hi') 



K 



Q(k-k')\/— + G(k!-k) 



(32) 



where k = k 2 , which is equivalent to expression (27) of 16]. 
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